Density-functional theory of the nonlinear optical susceptibility: Application to cubic semiconductors.

We present a general scheme for the computation of the time dependent (TD) quadratic susceptibility (χ(2)) of an extended insulator obtained by applying the ‘2n+ 1’ theorem to the action functional as defined in TD density functional theory. The resulting expression for χ(2) includes self-consistent local-field effects, and is a simple function of the linear response of the system. We compute the static χ(2) of nine III-V and five II-VI semiconductors using the local density approximation(LDA) obtaining good agreement with experiment. For GaP we also evaluate the TD χ(2) for second harmonic generation using TD-LDA. 42.65.Ky,71.10.+x,78.20.Wc Typeset using REVTEX 1 Nonlinear optics is a growing field of research which has applications in many technical areas such as optoelectronics, laser science, optical signal processing and optical computing [1]. In these fields the description of several physical phenomena, such as optical rectification, wave-mixing, Kerr effect or multi-photons absorbtion, relies on the knowledge of the nonlinear optical (NLO) susceptibilities. Moreover nonlinear spectroscopy is a powerful tool to analyze the structural and electronic properties of extended and low dimensional systems. In the present work we give a general scheme to compute from first principles the time dependent (TD) quadratic susceptibility (χ) of real materials within TD-density functional theory (DFT). Futhermore we show that the values of the static χ obtained in the local density approximation (LDA) are in good agreement with measured values for the cubic semiconductors. Our approach makes feasible the computation of χ in cells containing up to an hundred atoms, since it requires the same numerical effort as the computation of the total energy. This allows the evaluation of χ for systems of technological and scientific relevance which can not be handled by the traditional methods, such as surfaces or crystals of organic molecules. Nowadays many first-principle calculations for the ground state properties of materials are performed within DFT. Even in its simplest form, namely in the LDA for the exchange and correlation energy this scheme gives results which, in many cases, are in surprisingly good agreement with experiments. A rigorous extension of DFT to TD phenomena has been proposed in Ref.s [2,3]. Although the available approximations for the exchange and correlation energy are less accurate in the TD domain than in the static case, this scheme is sufficiently general to allow many possible improvements in the future. Therefore TD-DFT seems to be a promising framework for the study of the NLO susceptibilities. Standard quantum-mechanical perturbation theory can be used to compute the χ. The straightforward application of perturbation theory leads to an expression for χ, which diverges for an infinite solid in the static limit. However, for an insulator, these divergences have been shown to be apparent [4]. This kind of approach has been applied to compute the χ from first principles. The non self-consistent expression for χ reported in Ref. [4] has 2 been evaluated by Huang and Ching [5] using the DFT-LDA wavefunctions and eigenvalues. A fully self-consistent theory of the NLO susceptibility within DFT has been proposed in a series of papers by Levine and Allan [6]. Their method is feasible but algebraically very involved due to the necessity of dealing with the second order perturbation of the wavefunctions and with the apparent divergences. Their final expression is not easy to handle and its evaluation requires summations over the conduction band states, which are time consuming and difficult to converge. In a previous paper two of us [7] have shown that it is convenient to regard the static χ as a third order derivative of the total energy with respect to an uniform electric field. We pointed out that this derivative can be obtained by combining a Wannier representation of the electronic wavefunctions with the ‘2n + 1’ theorem of perturbation theory [8,9]. We also found an equivalent expression of the static χ in terms of Bloch wavefunctions. In the present letter we show that the method of Ref. [7] applies also to TD periodic perturbations and to the self-consistent TD-DFT functional. The TD χ can be regarded as a third order derivative of the total action. The stationary principle for the action functional [2,3], which replaces in the TD case the miminum principle for the energy functional, allows the use of the ‘2n + 1’ theorem. As in the static case the third order derivative depends only on the unperturbed wavefunctions and on their first order change due to the TD electric field. All the self-consistent contributions are included in the formalism in a simple way. The final expression avoids perturbation sums and does not present any apparent divergency. We apply our formalism to the computation of the the static χ of nine III-V and five II-VI cubic semiconductors within the LDA. For GaP we also evaluate the TD χ for second-harmonic generation (SHG) using TD-LDA [10]. In the Kohn and Sham (KS) formulation of DFT the ground state density n(r) of a system of N interacting electrons in an external potential Vext(r) is written in terms of N/2 single particle wave-functions {φ}. The set {φ}minimizes the KS energy functional E[{φ}] and the ground state energy is obtained as E = E[{φ}]. A formalism similar to that of the static case can be introduced also in the TD domain if one restricts to Hamiltonians periodic 3 in time and to the evolution of the system which is steady and has the same periodicity of the Hamiltonian [11]. In TD-DFT the TD steady density n(r, t) of a system of N interacting electrons in an external TD potential Vext(r, t), periodic in time with period T , is expressed in terms of aN/2 TD single particle wave-functions {ψ} [2,3]. The set {ψ}make stationary the KS action functional A[{ψ}], i.e. δA[{ψ}]/δ〈ψk(t)| = 0, (1) and the steady action is obtained as A = A[{ψ}]. The KS action functional A[{ψ}] is defined as (atomic units are used throughout) : A[{ψ}] = ∫ T 0 dt T 